library(haven)
setwd("~/Dropbox/FakeNews/final/analysis") # Change this with the path to the folder where you keep the data
rm(list=ls())

# Load the data previously saved by the Stata script file "fakeNews_script.do"
data <- read_dta("./data_fakenews_sim.dta")

#============================================================================================================
###### SIM1 - duplicate cases to make the interaction significant for Group 2: ANTI-GOV

output1 <- matrix(NA,300,100)

set.seed(1234)
for(i in 1:100){
  data1 <- data[data$typen==2,]
  for(j in 1:300){
    data1 <- data1[sample(nrow(data1)),]
    data1<- rbind(data1, data1[1:250,])
    lm1 <- lm(plaus ~ as.factor(fake)*polknow_irt*gov_eval + edu + age + male,data=data1)
    output1[j,i] <- coef(summary(lm1))[,"t value"][11]
  }
  cat("sim ",i," done \n")
}

final_output1 <- as.data.frame(cbind(seq(250,75000,250),0))
names(final_output1) <- c("morecases","avgt")

for(i in 1:300){
  final_output1$avgt[i] <- mean(output1[i,])
}

plot(final_output1$morecases,final_output1$avgt,main="Average t-score changes (MC)",
     xlab="Additional cases ", ylab="Average T-score", pch=1)
abline(1.96,0,col="red")

write.csv(final_output1,"./output_AntiGov.csv")

#============================================================================================================
###### SIM2 - duplicate cases to make the interaction significant for Group 3: PRO-GOV

output2 <- matrix(NA,300,100)

set.seed(1234)
for(i in 1:100){
    data1 <- data[data$typen==3,]
    for(j in 1:300){
        data1 <- data1[sample(nrow(data1)),]
        data1<- rbind(data1, data1[1:250,])
        lm1 <- lm(plaus ~ as.factor(fake)*polknow_irt*gov_eval + edu + age + male,data=data1)
        output2[j,i] <- coef(summary(lm1))[,"t value"][11]
    }
    cat("sim ",i," done \n")
}

final_output2 <- as.data.frame(cbind(seq(250,75000,250),0))
names(final_output2) <- c("morecases","avgt")

for(i in 1:300){
    final_output2$avgt[i] <- mean(output2[i,])
}

plot(final_output2$morecases,final_output2$avgt,main="Average t-score changes (MC)",
     xlab="Additional cases ", ylab="Average T-score", pch=1)
abline(1.96,0,col="red")

write.csv(final_output2,"./output_ProGov.csv")